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Abstract 



We consider a family of random matrix ensembles (RME) invariant 
under similarity transformations and described by the probability density 
P(H) = exp[— Try(H)]. Dyson's mean field theory (MFT) of the corre- 
sponding plasma model of eigenvalues is generalized to the case of weak con- 
fining potential, V{e) ~ ^ ln^{e). The eigenvalue statistics derived from MFT 
are shown to deviate substantially from the classical Wigner-Dyson statis- 
tics when A < 1. By performing systematic Monte Carlo simulations on the 
plasma model, we compute all the relevant statistical properties of the RME 
with weak confinement. For Ac ~ 0.4 the distribution function of the energy- 
level spacings (LSDF) of this RME coincides in a large energy window with 
the LSDF of the three dimensional Anderson model at the metal-insulator 
transition. For the same Ac, the variance of the number of levels, (n^) — (n)^, 
in an interval containing (n) levels on average, grows linearly with (n), and 
its slope is equal to 0.32 it 0.02, which is consistent with the value found for 
the Anderson model at the critical point. 

PACS numbers: 71.30.+h, 72.15.Rn, 05.60.+W 



1 



Typeset using REVT^ 



I. INTRODUCTION 



Random Matrix Theory (RMT) was introduced by WignerGi and DysonB to provide a 
statistical description of the quantized energy levels of heavy nuclei, and since then it has 
been applied to a great variety of complex systems, quantum and classical!. Gorkov and 
Eliashebergi suggested that the Wigner-Dyson (WD) statistics, derived from RMT, could 
be used to describe the energy levels of small metallic particles at low temperature, in 
connection with the study of their electromagnetic properties. Here a statistical description 
is made necessary by the presence of disorder and irregularities in the shape of the particles. 
For the case of disordered conductors one can resort to powerful field-theoretical techniques, 
which have allowed Efetovi and Altshuler and Shklovskiii to show analytically that the WD 
statistics are more than a simple phenomenological conjecture and describe exactly the local 
fluctuations of the energy levels in metals in a certain regime. 

The WD statistics are characterized by strong energy-level correlations, giving rise to 
the phenomenon of the level repulsion. These are the correlations that typically exist among 
the eigenvalues of a Gaussian Ensemble (GE) of random Hermitian matrices H, that is, an 
ensemble of matrices randomly distributed with probability density -P(H) oc exp[— TrH^]. 

The GE's (and therefore the WD statistics) do not bear any hint of the spatial dimension- 
ality (i of a physical system. Furthermore, they are, by definition, invariant under similarity 
transformations and thus there is no basis preference in them. This means that they can be 
applied only to particular regimes of a physical system where: 1) all the normalized linear 
combinations of the eigenstates have similar properties; 2) the dimensionality is, in some 
sense, irrelevant. 

For a disordered electronic system this is just the ergodic regime of the metallic state. 
In the metallic phase the eigenstates are extended structurless objects. If we further assume 
that all the relevant times are larger than the ergodic time r^) = L'^/D,{D being the diffusin 
coefficient), any diffusive particle can completely and homogeneously fill the total sample 
volume V = L"^ during its trajectory (ergodic regime) and thus does not feel the space 
dimensionality. Alternatively, the ergodic time defines a natural energy scale, Ec = h/ro 
known as Thouless energy. The ergodic regime, and therefore RMT, is valid within energy 
intervals e << E^, or within an interval containing a number of levels < E^/A ~ g, where 
A is the mean level spacing. The quantity g is known as the dimensionless conductance 
in units of e'^/h. When e > Ec the level statistics depends on the dimensionality and is 
different from WD distribution (Altshuler-Shklovskii regime). 

The non-ergodic regime is never reachable in the metallic phase in the thermodynamic 
limit, in any energy interval containing a large but finite number of levels, since the dimen- 
sionless conductance g, diverges in the limit L —>■ oo. Therefore the WD statistics describe 
exactly the energy- level correlations of the metallic state in the L — oo limit, that exists 
for d > 2 at relatively small disorder. 

For d = 1, 2 no metallic state exists in the thermodynamic limit, if any disorder is present, 
the system being always an insulator. For d = 3, with disorder increasing, the system goes 
through the Anderson transitionll to the insulating state, where all states are localized.! The 
level statistics in these situations obviously cannot be described by the a- dimensional, clas- 
sical RMT of the GE. In particular, such simple U(N) invariant random matrix ensembles 
(RME) cannot be an appropriate description of an Hamiltonian matrix whose eigenvec- 
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tors undergo the phenomenon of locahzation, since one can construct extended states by 
a hnear combination of locahzed states. In contrast, the proper -P(H) distribution should 
contain eigenvector-eigenvalues correlations or a basis preference that exclude those unitary 
transformations which would lead to the formation of such extended states. 

Random banded matricesi0 are an example of non invariant RME and perhaps more 
realistic for describing Hamiltonians of quasi-onedimensional disordered electronic systemsil 
and other quantum chaotic systems.0 Their statistical properties exhibit a crossover from 
the WD to the Poisson statistics as a function of the parameter h/N, [b is the band width) 
which is similar to what happens in quasi-onedimensional systems upon decreasing the ratio 
C,/L, ^ being the localization length. 

While the RMT description of the quasi-onedimensional disordered systems can be pro- 
vided by random banded matrices, an analogous description of the energy statistics near the 
critical point for c? > 2 is still missing. Thisnroblem has recently become a very outstanding 
one, after intense study, both numerica and analytical,liHlEy has shown that the spectral 
correlations at the metal-insulator transition are also universal and very different form the 
WD and the Poisson statistics. 

There is not yet full consensus on the exact nature of the critical energy-level statistics. 
The main finding of the analytical treatment based on a scaling analysisEl'ls is the existence 
of a power-law decay of the two-level correlation function with a non trivial critical expo- 
nent. All numerical simulations show that at the critical point the statistical fluctuations of 
the energy levels are scale invariant. However while some authors§'0@ claim a good agree- 
ment of their numerical results with the analytical predictions, othersS'EUll suggest that 
the critical statistics constitute more simply a "hybrid" between the WD and the Poisson 
distributions without any nontrivial critical exponent. 

The search of a RMT description of the critical statistics has prompted the investigation 
of physically motivated RME which exhibit non trivial deviations from the WD statistics. 
One important generalization has been obtained in Refs. [ ^ starting from the GE and 
introducing a symmetry breaking term of the form 

P(H) OC e-^'-H' e-/.^^^T.([A,H][A,H]t)_ ^^ -^^ 

The /i-dependent term breaks the U (N) invariance and tends to align H with a symmetry 



breaking unitary matrix A thus setting the basis preference. It was shown in Ref. [ ^ that 
even after averaging over A the implicit presence of the symmetry breaking term causes a 
dramatic change in the level correlations of the resulting ensemble, which exhibit a cross-over 
from the WD to the Poisson statistics. 

On the other hand, a lot of workU has been devoted to the analysis of U{N) invariant 
generalizations of the GE, of the kind 

P(H)oce-T'-^(«) (1.2) 

Such ensembles arise from a Global Maximum Entropy Ansatz of RMT in which an 
information entropy is maximized by the distribution. §1 The function V(H) acts like a 
generalized Lagrange multiplier and it is determined, e.g., by requiring that the density of 
eigenvalues is some given function /3(e), taken directly from the microscopic system being 
investigated: 
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(Tr[5(e-H)])=p(6) (1.3) 

For a long time it was believed that the the local statistical properties of the eigen- 
values of these ensembles are completely independent of V and identical to those of the 
GE. This hypothesis, so far supported only by numerical evidence, has been proved more 
rigorously very recentlyBS a large class of functions V. On the other hand, it has 
also been demonstrated&^ that there exists another class of functions V, known as weak 
confining potentials, for which the eigenvalue statistics display very strong deviations from 
the universal WD behavior. Such potentials are characterized by a very slow asymptotic 
growth: 

\/(e) ~ Aln^ |e| , |e| ^ oo . (1.4) 

It is important to emphasize that this asymptotic behavior of V has been inferred from 
numerical studies of Random Transfer Matrix models for disordered conductors, through 
the Maximum Entropy Ansatz. In that case one considers the eigenvalues of some 

combination of the transfer matrix, which are directly related to the conductance and become 
larger and larger, namelv less confined, when the disorder increases. 

It has been shownc3Ell that the local eigenvalue fluctuations of the RME with confining 
potential ( |1.4| ), exibit a cross-over from the WD statistics to a more Poisson-like behavior 
when the parameter A is decreased. Since it is a common belief that there is a connection 
between the statistics of eigenvalues and eigenstates of a RME, this breakdown of the WD 
universality seems to contradict the argument, (presented above for the GE), that a U{N) 
invariant ensemble cannot exhibit Poissonian statistics, typical of localized states. In Ref. [ 



3T| it has been suggested that Poissonian behavior in such RME is a remarkable phenomenon 
due to the spontaneous breakdown of the U{N) symmetry at the transition from a power-law 
potential, ^(e) ~ |e|°, to the logarithmic potential of Eq. (|1.4|) when a ^ 0. A crucial point 
in reaching this conclusion was the observation that the two-level correlation function of the 
invariant RME with weak confinement and the one of the RME with symmetry breaking 
are identical in a certain range of the parameters. 

Having established that these generalized RME with soft confinement belong to a new 
universality class, characterized by nonclassical correlations that interpolate from the WD 
to the Poisson statistics, there remains the important question of whether or not their 
properties are related to the critical energy level statistics of the Anderson model. 

In this paper we address this question through a careful study of the RME with dis- 



tribution given by Eq. [|1.3|| . Firstly we shall show that the local statistical properties of 



a generalized RME of this sort can be correctly and easily determined using an extension 
of Dyson's functional derivative formalism for the corresponding one dimensional plasma 
model for the eigenvalues. This formalism has been used successfully by Beenakkeiil to 
investigate the random-matrix theory of the transmission eigenvalues in disordered conduc- 
tors. In particular we prove that a modified mean field theory (MFT) of the plasma model 
is able to yield the nonuniversal asymptotic behavior of the two-level correlation function 
for weak confinement. 

Secondly, we perform extensive Monte Carlo simulations of the onedimensional Coulomb 
plasma of the eigenvalues and calculate all the relevant quantities that describe short-range 
and long-range statistical properties of the RME. The comparison with the MFT and other 
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analytical results shows that this Monte Carlo method provides very accurate answers for 
this problem and can be used to study more complicated RME where no analytical results 
are known. 

Finally we critically compare the level statistics of these RME with the the results re- 
cently obtained for the 3D Anderson model at the metal-insulator transition. Our analysis 
shows that, while the asymptotic correlations of the RME do not agree with the analytical 
result of Ref. [ |T^, two other statistical properties, namely the distribution function of the 



level spacings and level number variance, are remarkably close to those found numerically 
in exact diagonalizations of the Anderson model. 

The paper is organized in the following way. In Section we set up formalism and 
notation and review Dyson's derivation of the effective plasma model for the eigenvalue 
distribution. In Section [[TTTI we develop a MET of this model, generalized to the weak 
confining potentials. The MET study of the two-level correlation function is carried out in 
Section [[IV| . Section [|V| is devoted to the study of Monte Carlo simulations. In Section JVJ 
we discuss and compare the results of RMT with the the analytical and numerical results of 
the energy level statistics of the Anderson model at the metal-insulator transition. Summary 
and conclusions are presented in Section WW- 



II. EIGENVALUES STATISTICS AND EFFECTIVE ONE DIMENSIONAL 

PLASMA MODEL 

We consider an ensemble of random N x N matrices H. The matrices H are supposed 
to represent, for example, the Hamiltonian of a complex system, such as a quantum dis- 
ordered conductor. We take H to be either real symmetric, Hermitian, or quaternion-real 
self-duaM This choice defines three possible ensembles corresponding to three different phys- 
ical systems: 1) systems with time-reversal and rotational invariance; 2) systems with broken 
time-reversal symmetry; 3) systems with time-reversal symmetry but broken rotational in- 
variance. 

According to the maximum entropy principle mentioned in the previous section,!! we 
will assume that the probability distribution for the RME is defined by the density 

P(H) = Z-ie-T'-^(«) (2.1) 

where Z is a normalization constant. The volume element is d[Ii] = Yh>jdHij, for real 
symmetric matrices, with obvious generalizations for the other two cases.Q The probability 
density -P(H) d[Ii\ is evidently invariant under orthogonal, unitary or simplectic transfor- 
mations respectively, according to the three possible choices of H. The three ensembles are 
therefore specified and denoted after their internal symmetry: Orthogonal Ensemble (OE), 
Unitary Ensemble (UE) and Symplectic Ensemble (SE).I 

The invariance of -P(H) in Eq. implies that different matrices with the same eigen- 
values have the same probability to occur in the distribution. 

One can take advantage of this property and obtain the joint probability distribution 
(JPD), V{{ei}), for the eigenvalues e^, i = 1,2, . . . ,N of the matrices H. For this purpose 
it is necessary to express the various components of H in terms of the N eigenvalues and 
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other mutually independent variables pj, which together with form a complete set. The 
variables pj can be integrated out and the final result for JPD iJl 




nU^}) = CN,exp{-^PY.yi^)]]l (2-2) 



where /? = 1 for the OE, /3 = 2 for the UE and /? = 4 for the SE. Cn^ is such that V is 
normalized to unity. The universal Jastrow factor ni<fc — ^k\^ comes form the Jacobian 
of the variable transformation. It is universal in the sense that 1) it is always present in 
the JPD of the eigenvalues, whenever the the initial RME distribution probability is of the 
form Eq. (|2.1| ), 2) it is independent of the particular choice of l^(H) and depends only on 
the symmetry of the ensemble. 

The JPD of Eq. (|2.2|) characterizes all the statistical properties of the eigenvalues of 
an invariant RME, with the internal symmetry discussed above. It describes the so called 
energy-level statistics, if H is the Hamiltonian of the system. 

Following DysonE3, the JPD can be rewritten in the following form: 



V{{t,}) = Z^l exp[-/57^({e,})] (2.3a) 
'H{{en}) = - E In |e. - I + E ^(^0- (2.3b) 

i<j i 

The probability distribution ( p.2|) has the form of a Gibbs distribution for a classical, one 
dimensional system of "particles" e^, described by the "Hamiltonian" Ti. The symmetry 
parameter (3 plays the role of the equilibrium "temperature" . 

These fictitious "particles", namely the eigenvalues e^, interact among each other through 
a pairwise logarithmic repulsion at any energy scale. The external one-body potential V{e) 
keeps the system confined. V{e) is the only quantity of the RME that can be related to the 
microscopic parameters of the original physical model through its global statistical property, 
namely the density of eigenvalues, [see below]. The logarithmic repulsion does not depend 
on any microscopic detail of the real system, its origin being completely geometrical. From 
now on we will use freely this particle model analogy and call "particles" the eigenvalues of 
the RME. 

The local statistical fiuctuations of 1 ~ n << eigenvalues are conveniently described 
by the n-level correlation functions, defined as 

jyi /-oo /-oo 

gn(ei,...,en) = • / .../ VN{{ei))dei...den, (2.4) 

[JS — ny. J-oo J-oo 

By definition, Gn{^i, . . . , e„) is the probability of finding simultaneously any n particles 
at positions ei, e2, . . . , e„, the positions of the remaining N — n remaining unspecified. G„ 
are positive defined. In particular Gi(e) gives the density of particles at position e, and it 
will be denoted by 

p(e) = Gi(e) (2.5) 

It is convenient to introduce the n-level cluster functions or cumulants. The normalized 
rz-level cluster function is defined by in the usual way of statistical mechanics. The first two 
cumulants are 
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>^ ei, £2=1- . ^ . ^ , 2.7 
P(ei)p(e2) 

We shall always consider the case of large > > 1 . In this limit the normalized cluster 
functions are very useful since they tend to definite limits when the variables are written 
in the correct units. In taking the N ^ oo limit it is necessary to measure the particle 
positions in terms of the mean level spacing, A. If limAr^ooP(e) = Po = constant, A = p^^ 
and the dimensionless variables are simply 

Si = ei/A (2.8) 

On the other hand if limAr^ooP(e) 7^ constant, we need to consider a rescaling of the with 
the local density or the more complicated unfolding procedure [see Sec. [IV B| ]. 



In any case, the Si will form a statistical model for an infinite number of particles with 
mean spacing equal to unity. It is only when written in terms of these rescaled variables 
that the local statistical properties of the eigenvalues of different RME can be meaningfully 
compared. 

We now come to the discussion on the explicit form of the potential V{e) and how the 
statistics depend on it. The case of the Gaussian Ensembles (GE), with potential v{e) = 
is the only one for which exact solutions for the density, two-point correlation function and 
other statistical properties have been known exactly for a long time.i They are usually 
referred as Wigner-Dyson (WD) or classical statistics. One can show that in the large 

limit the particle density obeys Wigner's semicircle law, with a radius proportional to 
A^^/^. The two-point cluster function, in a region around the origin e = 0, is translationally 
invariant, ¥2(31, S2) = ^2('"), ^ = |si — S2I. For P = 2 (GUE) it has the famous form 



>2(r) 



sm vrr 



vrr 



(2.9) 



Similar expressions hold for the other two GE's.S 

At small separations the correlation function G2{r) = 1 — Y2{r) vanishes, due to the 
phenomenon of level repulsion brought about by the logarithmic interaction 

G'2(r)~r^, r«l (2.10) 

Although the global statistics (such as the density) of the energy spectra of real systems 
do not follow the semicircle law, the local statistics of the level correlations of many chaotic 
and complex systems are very well described by Eq. ( plol) . In particular Eq. (|2.9|) describes 
the correlations in small metallic samples at low temperature,@ as well as the correlations 
amone a large but finite number of energy levels of a metallic system in the thermodynamic 
limit 

Untill recently it was believed that the form of the confining potential could only affect 
the density of eigenstates but not their local statistics, which would therefore be universal 
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and equal to the WD statistics of the Gaussian Ensembles. Such universality has been 
indeed proved recently in a rigorous way, for a large class ofpotentials which confine the 
system strongly, that is when V{e) increases faster than |e|.E§EZl 

In what follows we will consider two kinds of non Guassian potentials, which confine the 
system weakly, with the following asymptotic behavior: 

1) Power-law PotentiafliUl 

A 

V^(e) = -|e|", 0<a<l, |e| ^ oo (2.11) 

2) Squared Logarithmic PotentialSH 

V^(e) = -In^ |e|, |e| ^ oo (2.12) 

It is a legitimate mathematical interest to investigate these cases and see if the WD 
universality is preserved. However these RME have also a physical intererst, since they have 
been suggested by random transfer matrix models of disordered conductors. 



III. MEAN FIELD THEORY (MFT) 

We now consider a mean field theory (MFT) analysis of the classical one dimensional 
plasma model of eigenvalues V{{ei}), into which the original probability density of RME has 
been mapped. We first define a continuous limit of this model,i@ valid in the asymptotic 
limit of large N, by introducing the particle density p(e) = this limit we will 

assume that the Coulomb gas is a classical fluid with a continuous and smooth macroscopic 
density. By substituting this definition of p(e) in Eq. ( p.3b| ) the Hamiltonian 7i{e„}) becomes 
an Energy functional, H[p{e)], of p(e) 

^ r+oo r+oo /■+00 

n[p{e)] = -- de de'p{e)p{e') In |6 -e'\+ dep{e)V{e) . (3.1) 

When N is large but finite, the assumption of a smooth density is only an approxima- 
tion. Consequently the first term of R.H.S. of Eq. (|3.1|) does not reproduce exactly the 
corresponding term of Eq. ( p.3b|) , because it neglects the two- level correlations, that is, it 
allows the presence of the "charges" p{e)de and p{e')de' at separations e — e' — 0. Since the 
interaction In |e — e'| is singular, this approximation has an effect, albeit small because N is 
large. One can compute the correction in a 1/A^ expansion to this term and see that it is of 
the form@ 

1 p+oo 

mP{^)] = -7; dep{e)Hp{e)) (3.2) 

Z J —oo 

It is convenient to introduce a grand canonical potential VL[p\ = 7Y[p] — pAf[p], where 
JV is the particle number functional A/'[p] = / p{e)de and p is the chemical potential. The 
average density of particle (p(e)) can be expressed in terms of the functional integral 

(p(e)) = Z-^ J p(e)e-^^[^]pp ] Z = j e-^^'^P^Vp . (3.3) 
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The one-body potential V{e) acts as a source term for the field p(e) and (p(e) can be 
expressed as 



(p(e) = ipZ)-'5Z/SVie) . (3.4) 

Up to now we have only assumed the existence of a smooth particle density (p(e), neces- 
sary in taking the continuous limit of the Coulomb plasma, which presumably is a good ap- 
proximation when N is large. The MFT, based on the continuous approximation, amounts 
to neglecting any entropy fiuctuations about the average and using the saddle-point ap- 
proximation in the integral equation for (p(e)). The MFT equation obeyed by (p(e))MF is 
thereforeii 

' rfe'(p(e'))MF In |e - e'\ = V{e) - p, (3.5) 

where the "Lagrange multiplier" p is to be found from the normalization condition 
I{p{^))uF de = N. The MFT interpretation of this equation is very natural, since it rep- 
resents the condition of mechanical equilibrium for the charge density (p(e))MF subject to 
the external potential V{e). Such a MFT approximation completely disregards the entropy 
part S[p] = — / de{p{e)) ln(p(e)) in the Free-Energy functional, jF[p] = H[p] — TS[p], and is 
exactly applicable only for 13 = oo. However, the long-range nature of the pair- wise interac- 
tion in Eq. ( p.l| ) makes the MFT approximation valid in the bulk of the spectrum even at 
finite f3, when N ^ oo. Indeed Dyson@ has calculated the first correction to this equation 
in a expansion and shown that the more accurate equation for (p(e)) reads 

/ de'{p{e')) In |e " ^'1 + ^ Hpie)) = V{e) - p, (3.6) 



The second term on the left-hand side of Eq. ( p.6|) is the sum of two parts. The term 
proportional to — T = — comes clearly from the entropy contribution to the Free Energy. 
The part independent of /3 is generated in the passage to the continuous limit as discussed 
above [see Eq. ( |3.6| )]. The combination of these two terms acts as an entropy contribution 
multiplied by an effective temperature, T*(/3) 

T*{(3) = f3-' - 1/2 , (3.7) 

which vanishes for (3 = 2. We can write the Free Energy, at equilibrium, of the system in 
the following form 

F = -l, 

2 J- 



rfe'(p(e))(p(e'))ln|e-e'|+ / de{p{e))V{e) 

r+oo 

T\(3) rfe[-(p(e))ln(p(e))]. (3.8) 



For the class of confining potentials of Eq. (|2.11|) and Eq. ( p.l2| ), one can show@ that 



the relative contribution of the correction to the MFT equation is of order A^~^lnA^, and 
therefore vanishes in the thermodynamic limit. For the time being we will neglect this 
correction and concentrate on the MFT Eq. (p.5|). 
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A. Solution of the eigenvalue density 



We now present the solution of the integral equation ( ^.51 ) for the particle density 
p(e)MF = (p(e))MF- The solution of this equation, confined to the region —D < e < D, 
can be found using the Cauchy method^ and is given by: 

where e+ = e + iO and the band-edge D is to be found from the normalization condition. 
We will discuss separately the two cases of a power-law and logarithmic potential. 

1. Power-law confinement, V{e) = ^|e|" 
The integral equation for (p(e))MF becomes 



= ^yi^Re^°_^_JL^, (3.10) 

We need to distinguish further between the cases a > 1 and a < 1. For a > 1, the integral 
in Eq. ( |3.1CI|) is divergent in the thermodynamic limit N ^ oo, when the band edge D 
also goes to infinity. The main contribution to the integral is made by the region ^ ~ D 
and therefore for any fixed e << D, one can neglect the e dependence in the integrand of 
Eq. (|3.10|) . Then the mean level density tends to a constant, p A^i~i/". Thus we reach the 



important conclusion that, for a > 1, there exists translational invariance in the e space in 
the N ^ oo limit, exactly as in the case of the Gaussian ensemble. Therefore for a > 1 we 
have a condition of strong confinement, and the corresponding local statistics belong to the 
WD universality class in agreement with Refs. [ On the other hand, for a < 1, the 

integral in Eq. ( ^.lU ) is convergent even in the limit D ^ oo. For finite ( and therefore 



finite D), the integral can be calculated by transforming the original contour of integration 
into a sum of two pieces, F = Fi U F2, where Fi is the negative imaginary axis, and F2 
is a part of the positive real axis, [D : 00]. Since we are interested only in the real part 
contribution, the integral over F2 does not contribute and we end up with 

-W-^-m^fTO?^ ,3.n) 

^"'-C.J!^r^Fa^-±^;l,-A, ,3.12) 



TiD ' \2J k^-" V2 2 2 



= ^, . =\ + 0{a') , (3.13) 

1 [a)L [1 + a) 2 

where z = e/D and F{a, b; c; x) is a hypergeometric function. The band edge is found from 
the normalization condition to be: 
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D 



l/a 



(3.14) 



The limiting function Pmf(^) — liii^A^^oo Pmf{^) is immediately obtained from Eq. (|3.12|) by 
taking the limit z = e/ D ^ 0: 



PMF(e) 



Aa /an 
tan 

277 



|1-Q 



(3.15) 



Therefore, for a < 1 the mean density in the thermodynamic limit, PmfI^)) does not scale 
as a power of N and, especially, is not translational invariant, p^p is divergent at e = 0. 
However, as we will show in Sec. 0, this is an artifact of the MFT solution and the exact 
density is sharply peaked but finite at the center. We can say that the system undergoes a 
"phase transition" at a = 1 and the symmetry which is broken is the translational invariance 
of the problem in the limit N ^ oo. 



2. Squared logarithmic confinement, V{e) = ^ln^(|e|) 

In order to solve the integral equation for the MFT density, it is convenient to represent 
this potential as a limit of combinations of the power-law potentiahil 



^(e) =ln2|e| = lim[a-2(|e|" - 1) 

a— >0 



(3.16) 



The MFT density in this case is 



pMF[e) = — Vl - z"^ 

TlJJ 



i^(|,| + «;|;l-r 



A-2a 



V2' 2 ' 2' 
"/2 ZT^ 



(3.17) 



where D satisfies the equation 



a 



C„ - 



2a 



r 1^ 



(3.18) 



The expression for the MFT density for N finite takes a simpler expression after the limit 
a — > is taken 



pMF(e) = A- 



arcsin \ 1 



The band edge is now an exponential function of A^: 

D = 2exp{N/A) 



(3.19) 



(3.20) 



The singularity at e = is again an artifact of the MFT solution. On the other hand we will 
see that in the bulk of the spectrum the mean density is accurately described by Eq. ( ^.19[ ). 
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IV. THE TWO-LEVEL CORRELATION FUNCTION 



In this section we study the two-level correlation function using the MFT developed in 
the previous section. Some of these results have been already derived for /3 = 2, using 
the different method of the orthogonal polynomials. MFT will allow to generalize these 
results to any value of (3. 



A. MFT Integral equation for the two-level correlation function 

Within the continuos formalism of the previous section, let us consider the following 
definition of the density- density correlation function R2{e, e') 

Using the definition of p{e) one sees that i?2(e, e') differs from the two- level c/uster function 
¥2(6, e') defined in Eq. ( p.9|) of Section I by the singular self-correlation 

>^2(6,6') = ^5(6-6') -i?2(e,e') (4.2) 

The correlation function -R2(e, e') can be easily expressed in terms of a functional derivative 
of (p(eVL with respect to V{e) in Eq. (|3.3|) . By using the relation 6Q/6V{e) = p(e) one 
obtains!^ 

Within the MFT it is possible to write down an integral equation for the two-particle 
correlation function.^ By taking the functional derivative 5/5V{e) in Eq. ( p.5| ) and using 
Eq. (^4.3|) one obtains 

J de"p{e")MF p(e)MF R2{e, e") In |e' - e"| = -p-'5ie - e') + p-'Sfx/5V{e), (4.4) 

This important equation is sometimes used to claim and justify, within the approximation 
of MFT, the universality of the correlations in RMT, or rather their independence of the 
potential ^(e). The argument is usually the following. In the large N limit one implicitly 
assumes that, at least in the region of interest, the average density scales like and goes to a 
constant, (p(e)) = po- Iii this case, since the two body-interaction is translational invariant, 
the two-particle correlation function must be translational invariant as well. In particular, 
the variational derivative of the chemical potential, which is the only term that depends 
explicitly on the confining potential and is not translationally invariant, must vanished If we 
now introduce dimensionless variables s = e/A, rescaled by the mean level spacing A = p~^, 
one obtains an equation completely independent of ^(e). 

/3 J dsR2{s - s") In \s' - s"\ = -S{s - s') (4.5) 
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From this one concludes that the hmiting correlation function -R2(e — e') cannot depend on 
the choice of V{e) and it is universal. Its asymptotic behavior is easily found by solving 
Eq. ( |4.5| ) by Fourier transformation. The result is of course the universal WD behavior of 
the Gaussian Ensembles. 



This reasoning is crucially based on the assumption that the average density tends to a 
constant in the thermodynamic limit so that it can be simply rescaled away. But we have 
seen that there are cases in which, if the potential is weak enough, the density is not at 
all a constant in the N oo limit and in fact it decreases steeply with e. The system is 
not translationally invariant and we cannot carry out the same simple rescaling as before. 
Thus the proof of the universality of the correlations is not applicable, and one may expect 
R2{e,e') to be different from the WD form. It is often said that in this case MET breaks 
down and its equations become invalid. This statement is not completely correct. The MET 
still gives reasonable and in fact accurate results for the particle density in the bulk of the 
spectrum (as the comparison with numerical simulations will show). Therefore this theory 
should work also for the correlations, when properly applied. Of course we can no longer 
disregard the presence of a nonconstant density, which is responsible of the deviations from 
universality. Below we suggest a simple way out of this problem through the notion of 
spectrum unfolding. 

B. Spectrum Unfolding. Modified MFT equation for weak confinement. 

On comparing the correlation functions of different ensembles it is necessary to choose 
energy units such that the mean level spacing is equal to one. This is trivial when the 
average density is a constant p, since in that case A = and the convenient variable is 
s = e/A = ep. The presence of a non uniform density in the thermodynamic limit makes 
this linear rescaling impossible, even locally if the density is a rapidly varying function of 
the energy, as in the case of the logarithmic confinement. This is a common situation in 
the study of complex energy spectra, where one needs to subtract off the unwanted effects 
of a non constant average density, in order to analyze the fluctuations around the average 
density itself. This difficulty is handled by the so called "unfolding procedure" , namely the 
introduction of a new variable s = s(e), in terms of which the density becomes a constant. 
The variable that serves this purpose is the integrated density of states 




(4.6) 




(4.7) 



By particle conservation 



(p(e)) de = (p(s)) ds , 



(4.8) 



where (p(s)) is the average density in the new variable s. Erom Eqs. ( [4.7| ) and ( |4.8| ) it follows 
that 
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(pis)) = 1, {s) = l 



(4.9) 



The definition of the two particle correlation function in the unfolding variables becomes 



where is the inverse function of s(e). Using these definitions in Eq. ( [4.4|) the integral 
equation for R2{s, s') i 

J ds" R2{s,s")\n\es" -es'\ = -(3-^6{s - s) , (4.11) 

The MFT equation for {p{s)) reads 

/ ds'{p{s'))\n\es-e,,\ = V{es)- p. (4.12) 

Joo 

Eqs. (|4.11| ) and (|4.12|) together mean that l3R2{s, s') is the solving kernel oS 

/ ds ip{s')\n\es — es'\ = (p{s) + const , (4-13) 

J oo 

that is 

^P{s) = / ds (3R2{s, s')ip{s') . (4.14) 



The additive constant in Eq. ([4.13|) has to be chosen such that dsipi^s) = since the 
variations in (p(s)) must occur at constant 

Eq. ( [4.11| ) is the integral equation for the correlation function i?2(s, s') in case of weak 
confinement. From this equation we can see that, in contrast with Eq. ( [4.4|) valid for strong 
confinement, now R2{s,s') depends implicitly on the original confining potential through 
the unfolding variable s(e). Thus it will not, in general, be universal. 



C. Solution of the MFT integral equation 



The solving kernel of the integral equation ( |4.13| ) must be found separately for each 



unfolding function e(s), which in turn depends on the confining potential. We are interested 
in the correlations in the bulk of the spectrum, where MFT works and thus we can use the 
MFT expression for the density in order to obtain the unfolding function. 

We first consider the case of the power-law potential with < a < 1, which is the 
simplest case of weak confinement where, in principle, non classical correlations can arise. 
Using Eqs. ( p. 15] ) and ([4.7]), the unfolding function in the bulk of the spectrum, in the large 
N limit is 

sie) = Asgn(e) rfe'-^ = Me^ sgn{e) , A = - tan (^— j , (4.15) 

or 
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e(s) = A-^/"|s|^/"sgn(s) (4.16) 
In order to simplify the derivation, we choose a = , k = 1,2, . . .. In this case 

/ S \ (2fc+l) 

<s) = [j) . (4.17) 

The integral equation Eq. (|4.13|) is now 

/ ds^is') In - s'^^'^+i)] = ^{s) + const . (4.18) 

J oo 

We change variables ^(^fc+i) _ ^ ^ _ j^jk+r and we define 

= (sjx^ , (4.19) 

^{x) = ^{s,) , (4.20) 

1 -2fc 



in terms of which Eqs. ( |4.13| ) and ( [4.14] ) become 

dx'tl!{x')ln\x-x'\ = <^{x) , (4.22) 



ip{x)= dx' (3R2{x,x') ,^{x') (4.23) 



and the condition on %p{x) is now J^^dxip{x) = 0. This is exactly the familiar case of 
pure logarithmic interaction with strong confinement considered above [see Eq. ( [4.5| )]. The 
solution for the solving kernel R2{x,x') is the WD universal correlation given in Eq. ( [4.61) . 
Using the relation between R2{x,x') and R2{x,x') and going back to the original variables, 
s and s\ we obtain 



This correlation function is not translational invariant. However we restrict ourselves to the 
case in which both s and s' are in the bulk and As = |s — s'| << s. In this limit we can 
expand the previous expression in power of As/s. The correlation becomes translationally 
invariant but what is even more important is that we get back exactly Eq. ( [4. 61) , namely the 
universal WD correlation function. Thus we conclude, that at the MET level, the power-law 
ensemble has completely classical statistics, even when the translational invariance is broken 
in the thermodynamic limit and the unfolding procedure is necessary. 
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We now consider the case of the double logarithmic potential. The unfolding function in 
the bulk of the spectrum, |e| >> 1 is now 

A /"l^l ,1 A 

s(e) = 77 sgn(e) / de- = In |e| sgn(e) , (4.25) 
/ Jo e Z 

or 

e(s) = e2|^l/^sgn(s) ^ 2 sinh (^^j . (4.26) 

We plug this function into Eq. ( |4.13| ) and we perform the change of variable, sinh (^^^ = 
X, with the corresponding redefinition of tp, ip and i?2- Following the same procedure used 
for the power-law potential, we arrive at 



cosh (f ) cosh (^^ 
''^ ^ sinh^ (^) cosh^ (4^) ■ ^^-^^^ 



Again the correlation function is not translational invariant. If s and s' are both in the bulk, 
|s'| , |s| >> 1 and have the same sign, we obtain the translational invariant expression 



1 

Ti'A^P sinh2"("^) 



R2in){s,s') = -—^ -TT77Z7Y ^ (4-28) 



where the subscript (n) stands for "normal" part. If |s — s'| > A, the argument of sinh cannot 
be expanded and A does not scale away. Thus the correlations for the double logarithmic 
potential are no longer universal and in contrast with the power-law behavior of the WD 
class, they decrease exponentially in agreement with the exact solution for = 2 by Muttalib 
et alM 

However there is one more great surprise, pointed out in Ref. [ |31| for f3 = 2, which is 
the re-appearance of strong correlations at s' ^ —s. In fact when s and s' are in the bulk 
but have different sign, s' = —s + As we obtain 



where the subscript (a) stands for "anomalous" part. This anomalous part of the correlation 
function breaks dramatically the translational invariance. Its remarkable property is a nar- 
row correlation hole at s' ^ —s with a depth, controlled by A, that does not decrease when 
|s — s'l ~2|s| — >cxD . Notice also that the two regions where R2(n){s, s') and R2(a){s, s') are 
nonzero are separated by a very large distance when s and s' are in the bulk of the spectrum 
and is large. 

In Ref. [ a simplified application of the method of the orthogonal polynomials, valid 
for q = e~'^ << 1, was used to derive the two-level cluster function for the case P = 2, 
yielding the result: 



1 ^ cosh (^] cosh f ^ , 

y.M ^ ^ Is...,. - J^J^. . (4.30) 
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The normal part of l2(s, s') 



/ /N 1 sin7r(s — s') , 

^^'»''^-''' - ^ [-fr(iiiy] • f^-^^' 

is identical to the exact solution of Muttalib et a/S for the same small values of q and it 
also identical to the exact solution of Moshe et a/.i for the case of a RME with a symmetry 
breaking term. 

Here, using the MFT theory, we have generalized this result to any /3. As in all other MFT 
calculations, only the asymptotic form of the R2 is obtainable within the MFT treatment. 
In particular the oscillatory function which vanishes as |s — s'| ^0, and thus gives rise to a 
residual level repulsion, is totally out of reach. In its place we have 1/(3 which is the average 
of the oscillations. Thus MFT and the method of orthogonal polynomilas (which is exact) 
are completely consistent even in the case of weak confinement where deviations from WD 
occur. 

The appearance of anomalous correlations at s' ~ —s is the result of the system trying 
to develop Poissonian-like correlations at s' ~ s, while, at the same time, complying with 
the U{N) invariance which forces a normalization sum rule on R2 to be satisfied even in the 
large limit. Ill The sum rule reads 



R2{s,s')ds' = , (4.32) 

J —00 

or, in terms of ¥2(3, s'), 

roo 

Y2{s,s')ds' = 1 , (4.33) 



00 



and it is must be satisfied in the case of a U{N) invariant RME, because of the long-range 
nature of the universal logarithmic interaction always present in such ensembles.Ej The 
normal part alone of ¥2(3, s'), given in Eq. ( [4. 311 ), does not satisfy Eq. ( |4.32| ). But the sum 
rule deficiency 

roo 

V = l- Y2^,,){s)ds, (4.34) 

is taken care of by the anomalous correlational: 

00 

Y2i,)is) ds = ri . (4.35) 

o 

However, no matter how this is realized in practice, what is important is firstly the fact 
that the U{N) invariant RME with soft confinement manages to develop the exponentially 
decaying two-level correlation function of Eq. (|4.31| ) in an infinitely large energy region in 
the bulk of the spectrum, where the anomalous correlations are irrelevant. Secondly, the 
normal part of correlation function is exactly equal to the expression obtained for the RME 
( p..l| ), where the symmetry is explicitly broken. This occurrence has been interpreted in 
Ref. [ as a signal of the spontaneous breakdown of the U{N) invariance in case of soft 
confinement, with the parameter 7] playing the role of the order parameter. 
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The physical interpretation of the anomalous correlations at s' 



-s in terms of 



Coulomb plasma is very simple if one looks at how the logarithmic pair-wise interaction 
is transformed in the unfolded variables. 



f{s,s') = -ln|sinh(2s/A) -sinh(2s/A)| . 



(4.36) 



In the unfolded variables the interaction is no longer translational invariant. We can try to 
rewrite it in a form that looks more translationally invariant in the following way: 




cosh 



's + s'^ 



A 



In cosh 



A 



(4.37) 



We see that the interaction splits into two terms: the first is an ordinary, translationally 
invariant, repulsion between two particles located at s' and s; the second, however, represents 
the interaction between a particle at position s' with the image charge at —s of a particle 
at position s. In other words, a particle at position s will repel particles aroung it. But its 
image with respect to the origin will also repel particles around the the position —s. Notice 
that the image charge term depends on A - it increases when A decreases- but it is less 
singular that the direct term because of the functional dependence on cosh(s + s'). 

So far we have always restricted ourselves to bulk properties of the correlations. As 
we will show in more detail in the next section, there are some interesting effects in the 
correlation function in the center of the spectrum. There, even the power-law potential for 
a < 1 displays deviations from the universal WD behavior. 



V. MONTE CARLO SIMULATIONS 

The statistical properties of the onedimensional classical system@ in thermodynamical 
equilibrium, whose probability distribution is given by Eqs. ( p.3aj) -( p73E| ), can be conve- 
niently studied by carrying out Monte Carlo (MC) simulations. 

The MC method is useful because it is not restricted to a particular value of as 
for the method of orthogonal polynomials, and it is very accurate in all the regions of 
the spectrum, in contrast with MFT which is good only in the bulk. It also allows the 
evaluation of important statistical quantities like the level spacing distribution function 
(LSDF) and the number variance in a straightforward way. Finally it can be used to study 
level statistics where the "particle" interaction is more general than the simple logarithmic 
interaction considered in RMT. These more complicated interactions have been shown to 
play an important role in some disordered systems.0 

The nature of MC simulations is best illustrated with an example. Suppose we want to 
calculate the mean particle density p(e) as a function of the position e. According to Eq. ( p.5| ) 
we need to perform a multi-dimensional integral, which up to an overall normalization, is 
equivalent to a 'thermal' average over an ensemble of particle configurations. The MC 
method replaces this ensemble average with a 'time' average, but the time evolution is 
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determined by equations that are artificial and chosen for convenience. To calculate p(e) 
one partitions the real axis into bins with boundaries determined by 

e„ = nAe , n = ±1,2,3,... (5.1) 

where Ae is the with of the bins. At the end of each time step (to be defined below) we 
obtain an updated configuration {ej} ,i = 1, . . . ,N, and we add one to a bin if a particle 
in this configuration lies in it. As 'time' progresses, the number of particles in a bin will 
become proportional to the mean density density at the position where the bin is centered. 
For the evolution density we have taken a simple Metropolis algorithm, which works in the 
following way. At each time step or sweep, we scan through the particles and attempt to 
move each one. Actually in each sweep we pick times one particle at random in the 
system, so it is possible that in one particular sweep one particle is chosen more than once 
and another is not touched. The moving attempt involves picking at random any position 
between the particle the proceeds and the one that follows the particle that we are trying to 
move and taking this position as the new attempted position. The attempted move is chosen 
in this particular way simply to optimize the convergence rate of the algorithm. It has the 
important property that if we start with an order sequence of particles, ei < 62 < . . .e^, 
the sequence remains ordered in the time evolution. To decide whether or not to accept the 
move, we calculate the change AE to the system's energy that would occur if the particle 
were moved to the new position. If AE is negative the move is accepted, and the particle 
is given the new position. If AE is positive the move is accepted conditionally. One picks 
a random number between and 1 and accepts the move if this number is smaller than 
exp[—i3AE). Before measuring any quantity the system must reach equilibrium and this 
is obtained by a certain number of "warming-up" sweeps that reach a 'typically' sampled 
configuration. 

We have carried out simulations over systems with up to 200 particles. We noticed that 
the simulations are in all cases very stable even for smaller N, and therefore we have typically 
worked with systems of A^ = 100 particles. Equilibration is usually reached very fastly and 
we have typically used 10^ sweeps to warm up the system. The averaged are taken over 10^ 
sweeps and the statistics that we are able to obtain are usually excellent. 

To make sure that the method works and is able to give numerically accurate results 
we have first studied the three Gaussian Ensembles whose density, two-particle correlation 
function, spacing distribution and particle variance are exactly known.0. For these ensem- 
bles the MC simulations reproduce the known result very accurately. For example, in the 
calculation of -P(cr) the method is able to detect the small deviation that there is between 
the exact result and the Wigner surmise. 

We now proceed to discuss in detail the results for the different quantities of interest for 
the case of power-law and double logarithmic ensemble. 



A. The density of states 

The MC evaluation of the average density of states p(e) is carried out as explained in the 
example above. In Fig. |l]we plot this quantity for the logarithmic potential for A = 0.5 and 
/3 = 1,2. The case of weak power-law potential, with < a < 1 is qualitatively the sam 
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The agreement between the MC result and the MF expression pmf given in Eq. ( |3.19|) is very 
good except around the origin, where in contrast to MFT, the simulations give a sharply 
picked but finite density at e = 0. For (3 = 2 the MC result coincides with the value obtained 
by the method of orthogonal polynomialsEl The first hundred of these for the power-law 
potential can be easily generated numericallyl^i and the density at e = obtained from them 
converges very fast. This fast N independence of the center of the spectrum is also seen very 
well with the simulations and it is an important property of the particle density for weak 
confinement. We refer to it as the "incompressibility" of the core of the particle-density 
distribution. In contrast to the case of strong confinement, V{e) ~ |e|",a > 1, where the 
density at the center of the spectrum scales with A^"~^, for a < 1 the confining potential 
is too weak to "compress" the particle in the core region near the origin. After the initial 
formation of the sharp but finite peak at e = (which happens for N « 100), on adding 
more particles to the system, these always go to the ends of the distribution instead of 
spreading homogeneously throughout the spectrum. The particle density in the core of the 
spectrum is almost independent of but depends on the inverse temperature j3. In the 
bulk of the spectrum the density decreases like l/e^~° in the N oo limit. Therefore 
translational invariance is broken in the thermodynamic limit, in agreement with the MFT. 

B. The two-level correlation function 

The MC evaluation of the two-point correlation function R2{s., s'), as for any other corre- 
lation function, faces the complication of the breakdown of translational invariance. There- 
fore we need to carry out a numerical unfolding of the spectrum in order to compare with the 
classical statistics. We have considered three different unfolding procedures that can be used 
in different circumstances. In the simplest case we are interested in the bulk correlations. 
Therefore the simplest unfolding scheme consists in carrying out, for each MC configuration 
generated at "time" t, {ciY ,i = 1, . . . ,N, the mapping 

UiY {siY , Si = sgn(e) / de p(e)MF • (5.2) 

Jo 

The unfolded configurations {siY generated in this way are then used to measure the cor- 
relations, which will be automatically in the right units. 

The second method, which turns out to be particularly useful for the logarithmic po- 
tential, consists in performing the change of variable e —>■ s = s{e) directly in the joint 
probability density function 

^({e.}) ^ 'Pi{s^}) oc exp[-pni{s,})] , (5.3) 

ms..Y = -Elnks,: -e.J +EnesJ - ^Eexp[ln(l/p(6,J] . (5.4) 

Notice that the one-body confining potential is modified by a /^-dependent term coming 
from the Jacobian of the transformation. 
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Both these two unfolding schemes make use of the MFT particle density and therefore 
can only be used to study properties in the bulk. In order to study the correlation function 
around the origin we need a more precise expression for p(e). Therefore we carry out the 
following procedure. The unfolded two-level correlation function at s = can be expressed 

as 



i?2(0,s') =i?2(0,e,0 . (5.5) 

The function e(s) is now obtained by numerically inverting s(e) = /o^(p(e')) de' , with the 
density (p(e)) evaluated directly by the MC simulations. 

Once we have taken care of the unfolding, the evaluation the two-level correlation func- 
tion -R2(s, s') by MC is very simple: we fix a reference particle at s and then we compute 
the "conditional" particle density p{s')\s of all the remainingoarticles with respect to the 
reference point, using a partition in bins as explained before.E3 In studying bulk properties 
the reference particle can be let free to move, since in the unfolded coordinates the density 
is constant. In this case s in R2{s,s') must be interpreted as (s). Because of the particular 
way of choosing the trial moves adopted here, the standard deviation from (s) turns out to 
be small once the system has reached equilibrium. 

For the power-law potential the Monte Carlo simulations show that the two-level correla- 
tion function R2{s, s') is perfectly equal to the WD expression for any a when the reference 
particle is in the bulk of the spectrum, in agreement with the MFT and orthogonal polyno- 
mial results. However this universality is broken around the origin. In Fig. |^ we plot R2{0, s) 
for a = 0.2, and (3 = 2 (the 5 function is not included). For small s, -R2(0, s) does not follow 
the classical universal behavior but rather starts out like s^/"-. Thus we have a sort of 
"super Wigner" behavior, with stronger level repulsion at short separation the smaller is 
a. Since the cluster function must satisfy the normalization sum rule, this implies a faster 
decay of -R2(0, s) at large s. Using the method of orthogonal polynomials one can show that 
this decay goes like l/s'^^"'"^/"^0 

We now discuss the bulk properties of R2{s,s') for the logarithmic potential. As shown 
if Fig. 1^, the MC simulations confirms fully the surprising result of the bulk breakdown of 
the translational invariance in R2{s,s') and the appearance of the "ghost" correlation hole 
at s = —s'. The reference particle was let free to move in the positive part of the spectrum 
around (s) =24 (in unfolded coordinates). Besides the usual correlation hole around this 
position, another one appears symmetrically with respect to s = 0, as if there was an image 
of the reference particle, the "ghost", located around (s) = —24. The contribution of the 
anomalous part increases upon decreasing A. For values of the parameter A not too small, 
both parts of the cluster functions obtained numerically are in good agreement with the 
analytical results given in Eqs. ( 4.2^ ) and ( 4.29| ), as seen in Fig. ^ where the normal part 



of i?2 is plotted. However, as one can also see from the same figure, already for A = 0.2 
the MC result for the normal part starts to deviate from the analytical formula, which in 
fact becomes invalid for 74<l/7r^^0.1. In this regime the MC simulations show that 
the cluster function instead of decaying exponentially at very small s, starts out more and 
more flat and in the limit of very small A it converges to a box or "well" of width 1/2 and 
depth (-1). A similar behavior occurs for the anomalous part as well. Therefore, in the 
A —>■ limit, the cluster function is composed of two rectangular wells centered at s and —s. 
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This result can also be obtained from the exact solution by Muttalib et al.fB plotting their 
general expression for i?2, which remains valid in the regime of very small A. 

We conclude that the correlation function R{s, s') of the double logarithmic ensemble, 
displays a cross-over from WD toward a Poissonic behavior for intermediate A, i.e. I/tt^ < 
A < 1, but never really becomes exactly equal to 6{s — s') of the Poisson distribution. 

C. The Level Spacing Distribution Function (LSDF) 

Let us consider a sequence of successive levels ei < 62 < • • • and let Si, S2, ■ ■ ■ be their 
distance apart, Si = e^+i — e,. The average value of Si is the mean level spacing A. We 
suppose for the moment that the average density and therefore A = = const. We further 
define the relative spacings cxj = S^/A. The nearest neighbor LSDF, P{<7), is defined by the 
condition that P{a)ds is the probability that any cTj will have a value between a and a + da. 

If the energy levels are completely uncorrelated, one can immediately prove that the 
LSDF is the Poisson distribution, 

Pp((T) = exp(-a) . (5.6) 

In contrast, for a large class of chaotic or disorder systems where the energy levels are 
correlated, -P(ct) is very well described by the so called Wigner surmise,0 

PM = —exp{--a'). (5.7) 

The Wigner surmise vanishes at short separations, showing the phenomenon of level 
repulsion, typical, for example, for extended wave functions of a disordered conductor. The 
Poisson distribution, in contrast, allows level degeneracy, as in the case of an Anderson 
insulator, where the wave functions are localized and do not overlap. Notice that Pw{<^) falls 
down faster at large a than Pp{cr)- This is again due to level repulsion which, in a finite 
energy window, prevents the appearance of large energy gaps with no levels in them. 

For the Gaussian Ensembles one can derive exact expressions for -P(cr), which turn out 
to be very close but not identical to the Wigner surmise .i All the RME with logarithmic 
interaction and strong confinement belong to the GE (or WD) universality class and thus 
their P{(t) is also very close to the Wigner surmise. 

The analytical determination the LSDF is not straightforward. For (5 = 2 the function 
P((t) can be expressed in terms of a determinant of the two-level cluster function!, which 
in general must be evaluated numerically.il Alternatively one can use MFlliiSii but this 
method has not been extended yet to the potential (|2.12|) . On the other hand the LSDF is 
easily calculated by MC. In terms of the plasma model, -P(cr) is defined, once the system 
has been unfolded, as the probability density of finding the nearest adjacent particle at 
a distance s from a given reference particle. The LSDF obviously coincides with the the 
two-level correlation for very small a. 

From the results of the two-level correlation function we expect the LSDF for the power- 
law potential to be identical to the Wigner surmise in the bulk, with possible deviations at 
the origin. The MC simulations confirm fully these expectations. To calculate the unfolded 
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spacing around the origin, we fix a particle at e = and we perform tlie unfolding by 
computing 

(5.8) 

e=e{a) 

where the function e{a) is again obtained by numerically inverting a{e) = /o*^(p(e')) d-^', The 
result is shown in Fig. ^, where we plot the LSDF for aO.2 and (3 = 2. The classical spacing 
for the Gaussian Orthogonal Ensemble {a = 2) is also plotted. The figures clearly shows 
the deviations of -P(cr) from the GOE result. In particular for small a the LSDF does not 
follow the universal behavior o"^ of the Wigner surmise but starts out like cr^/". This is the 
same "super Wigner" behavior already found for the two-point correlation function. 

In Fig. D we show the bulk LSDF for the double logarithmic potential in the case of 
the Orthogonal Ensemble. We plot P(cr) for several values of the parameter A, together 
with the distribution of the Gaussian Orthogonal Ensemble and Poisson distribution for 
comparison. We see that, for those values of A < 1 for which the two-particle correlation 
function displayed a deviation from the classical GOE behavior, we have a corresponding 
deviation form the classical LSDF toward a more Poisson-like behavior. P(cr) still starts out 
linearly at small a, but the initial slope increases upon decreases A as a result of a smaller 
level repulsion. The peak of the distribution shifts from a ~ 1 to smaller values ~ A. For 
large separations the decay is also slower than the the GOE result. In fact, as we show in 
the inset, plotting In P{a) vs. a » 1, we can fit rather well the curves with straight lines 

\nP{a) -a{A,P)a , a»l, (5.9) 

where the constant a{A, (3) > 1 decreases with increasing A. Notice that all curve cross at 
the same point at a ~ 2. Similar features and deviations from the Gaussian Ensemble are 
obtained also from the other two symmetries, unitary and symplectic. 

The cross-over toward the Poisson distribution stops, however, at around A ^ 0.2. We 
have shown that for A < 0.2 the correlation function, instead of becoming closer and closer to 
a S function when A is further decreased, it turns around and for very small A it approaches 
instead a square well of width 1/2. Something similar happens for the LSDF. We can already 
see in Fig. |] that for A = 0.1 the initial slope of the distribution has stopped increasing, 
and the height of the peak is getting close to one. For yet smaller values of A (not shown 
in the picture), the initial slope starts decreasing and the LSDF, instead of approaching the 
Poisson distribution, will tend eventually to a single narrow peak of height >> 1, centered 
at cr ~ 1. 

D. The Number Variance 

So far we have considered the correlation functions that probe essentially the local 
fluctuations of a small number n ~ 1 of energy levels. We now turn to the variance 
var((n)) = {{n — (n))"^) = (n^) — (n)^ of the number of levels in an energy window that 
contains 1 <^ (n) <^ N on the average. The number variance is a statistical quantity that 
provides a quantitative measure of the long-range rigidity of the energy spectrum. 



P(cx) 



(p(e)) 
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For the Poisson distribution the levels are uncorrelated and there are large level-number 
fluctuations, leading to a linear variance 



varp((n)) = {n) . (5.10) 

On the other hand, the level correlations in the WD statistics makes the spectrum more 
rigid and the number variance grows only logarithmic 

var„((n)) = ln(n) + + 0(l/(n)) , (5.11) 

71 p 

where C/3 is a constant of order 1, which depends on the symmetry. 



The MC results for the power-law confinement show that Eq. ( |5.11| ) is perfectly satisfied, 
for every a, in the bulk, namely when the energy windows does not contain the origin with 
its nonuniversal correlations. On the other end when the energy windows are centered at 
the origin, the "super Wigner" correlations present at e = for a < 1, manifest themselves 
making the constant Cp in the number variance a dependent. As we show in Fig. |^, for 
f3 = 1, Ci decreases with a, when a < 1, because there is more level repulsion in the 
area of the origin and therefore more level rigidity in the overall spectrum. For a > 1 the 
universal value of the Guassian Ensemble Ci ~ 0.4420, is recovered. The a dependence of 
Cf3 is however the only deviation from universality in the number variance, the logarithmic 
dependence being unchanged for the power-law confinement. 

We now come to the case of the logarithmic potential. We have seen that the presence of 
the "ghost" correlations break translational invariance in the two-point correlation function. 



It was shown in Ref. [ ^ that, due to such a breakdown of translational invariance, 
the number variance for the double logarithmic confinement depends on the position of the 
energy window in a very essential way. If the energy window does not contain the origin, then 
the effect of the ghost peak is not felt and the system is Poisson-like, with translationally 
invariant correlations (in the energy range considered) given by Eq. ([4.28|) . In this case the 



number variance is also Poisson-like and increases linearly, as one expects in presence of 
exponentially decaying correlations. The coefficient of the linear term rj is less than 1, and 
it is given by Eq. ( [4.34|) . Here 77 is nonzero because the normal part alone of the correlation 
function fails to satisfy the normalization sum rule. It increases upon decreasing A because 
the "spectral weight" of the normal part of the cluster function ( [4.31|) decreases. However, if 



the energy windows are symmetric with respect to the origin, the ghost correlations become 
effective and their contribution allow the sum rule to be satisfied. Therefore the coefficient 
of linear term r] in the variance vanishes in this case. 

Indeed, the Monte Carlo simulations show a dramatic difference in the level number 
variance in these two cases. In Fig. ^ we show the number variance for an energy window 
centered at a point in the bulk, excluding the origin. The variance grows linearly and 
the coefficient of the linear term is in good agreement with Eq. ( [4.34|) and with the result 



obtained for the exactly soluble models by Blecken at al.E3 On the other hand. Fig. |^ shows 
the variance calculated for symmetric energy windows, containing the origin. The linear 
term is absent and the variance is constant for all integers (n) >> 1. Thus, despite the 
smaller level repulsion, the overall "level" rigidity is even higher that for the classical RMT. 
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VI. DISCUSSION AND CONNECTION WITH THE CRITICAL LEVEL 
STATISTICS OF THE ANDERSON MODEL 



Let us first summarize the main results of our analysis of the generalized RMT with soft 
confinement. We have seen that for the very weak logarithmic potential, ^(e) ~ Aln^ |e|, 
the local level fluctuations in the bulk of the spectrum display a cross-over from the WD to a 
more Poisson-like behavior, when the parameter A is decreased. In particular, the two-level 
correlation function, in the bulk of the spectrum far from the origin, decays exponentially at 
large distances. The spacing distribution function still vanishes like at short separation, 
but the initial slope is steeper, implying less level repulsion. The tail of the distribution 
decays like exp{—a{A,P)a) with a{A, (3) > 1, intermediate between the WD surmise and 
the Poisson function. The level number variance, when calculated within energy windows in 
the bulk of the spectrum that exclude the origin, is also Poissonian, increasing linearly with 
the the average number of levels (n). We saw, however, that the Poisson limit cannot be 
reached fully within this ensemble. The nonuniversal behavior of the RME with logarithmic 
confining potential has been attributed to a spontaneous breaking of the U{N) invariance. 
For steeper confining potentials, ^(e) = |e|", no deviation from the WD statistics occurs 
in the bulk of the spectrum, the only small deviations from universality occurring at the 
center .0 

The question that now we want to address is: does the non classical (namely non WD) 
behavior of the RME with double logarithmic confinement have anything to do with the 
universal energy-level statistics of the Anderson model at the critical point? A first point that 
needs addressing, before any comparison of the different statistical properties is attempted, is 
the way in which nonuniversality comes about in the invariant RMT with weak confinement. 
We have seen that the essential ingredient to obtain a deviation form the WD statistics is the 
strong energy dependence of the averaged level density: even in the N oo limit, (p(e)) is 
a rapidly varying function of e everywhere in the spectrum. In fact one cannot even define a 
constant local density, since the relative variation of (p(e)) over an energy range equal to the 
mean level spacing is of order 1. This is at odds with the well-known result for the density 
of states in the Anderson model: in that case p(e) is constant over a large energy region 
around the center, and moreover it is a noncritical quantity, that is, it does not exhibit any 
drastic change at the critical point. One can reply to this serious objection by reminding the 
similarly well-known fact that quite often complex systems with different alobal statistical 
properties (such the density of levels) have the same local level fluctuations. a and viceversa. 
The most famous example is the GE itself, whose semicircle law for the density of states is 
certainly not obeyed by any of the spectra of the heavy nuclei or other complex systems; 
yet its correlations are very universal and describe accurately the local statistical properties 
of these systems. Something of this sort might happen in our case. In this respect it is of 
great importance the existence of another ensemble, namely the broken symmetric model of 



Eq. (|1 . 1|) , where the level density is constant in the thermodynamic limit, and nevertheless 
the local level correlations are the same of the RME with weak confinement. Thus there 
exist at least two RME, having very different global statistics, whose local statistics belong 
to the same new universality class. We also must emphasize again that the asymptotic 
logarithmic behavior for the confining potential of the invariant RME has been suggested by 
studies on transfer matrix models of disordered conductors through the maximum entropy 
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principle.aH^ In the transfer matrix formalism one can express the conductance g in terms 
of the eigenvalues Xj of the matrix X = TT'^ + {TT^)~^ — 21, where T is the transfer matrix 
and / the unit matrix 

Localization appears in the presence of exponentially large eigenvalues Xi. Therefore a simple 
maximum entropy principle can provide, through the average density p{x), information 
about the localization of the wavefunctions, in contrast to what happens to the Hamiltonian 
matrix. Here we have assumed that an U{N) invariant RMT for the Hamiltonian matrix can 
be constructed from the confining potential derived from the corresponding transfer matrix. 
The hope is that such RMT will generate the correct local energy-level statistics, despite 
the aearge energy density itself will not be well represented.ii This procedure is probably 
too naive, but it is clearly the simplest and we will discuss its implications. 

The conjecture of the existence of universal statistical properties at the metal-insulator 
transition was put forward by Shkloskii et a/.Jill on the basis of numerical studies of the 
spacing distribution function (LSDF) which turned out to be scale invariant at the critical 
point. In a recent work Kravtsov et a/.El have carried out an analytical study of the critical 
statistics of the Anderson model. By using the analytical properties of the diffusion prop- 
agator and certain scaling relations valid at the mobility edge they have proved that the 
two-level correlation function has the following asymptotic behavior: 

R{s,s') = Cf3-^\s-s'\-^+^ , |s-s'|»l, (6.2a) 
7 = 1- (z/d)-i < 1 , (6.2b) 

where C is a positive constant while 7 is an universal critical exponent related to the critical 
exponent u of the correlation length ^. 

However the level number variance at the critical point contains two terms0 

var((n)) = ri{n) + b{ny , (6.3) 

where rj < 1 and b are some universal positive constants. The power-law term originates 
directly from the asymptotic power-law tail in the critical two-level correlator and thus 
reflects the critical dynamics. But there is also a linear term, which had been already 
predicted by Altshuler et a/.0. Formally the origin of this term is again due to the violation 
of the sum rule (f4.34|) , by the critical two-level correlation function. The physical meaning 



of this term is not yet understood. Its existence however implies that the dominant term 
in the variance at the critical point is still Poissonian, albei with coefficient less than 1. 
Notice that the knowledge of the two-level correlation function is not sufficient to develop 
a complete statistical description of the energy level at the transition. For example, the 
LSDF cannot be found without further statistical assumptions. By mapping the critical 
energy levels into a plasma model and assuming a the existence of a particular pair-wise 
interaction^, one can use the analytical result of Eq. ( |6.2a|) for the two-level correlation 



function to derive explicitly the effective repulsive interaction among the levels.0. Once the 
resulting interaction is known the asymptotic form of the LSDF can be evaluated, obtaining 
the result^ 
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P{a) ~ exp (-/i^a^-^) , (6.4) 

where h-y is a positive constant. Despite some numerical simulations seem to support this 
finding (but see below), this approach has the serious drawback that it gives rise only to 
the second term of Eq. (|6.3|) for the level number variance, the linear term being absent and 
totally unexplained. 

Following the work by Shklovskii et al, several other groups have studied numerically 
the energy level statistics at the critical point. In all cases the statistical fiuctuation prop- 
erty that is easiest to study numerically, namely the LSDF, shows scale invariance and a 
behavior intermediate between the WD surmise and the Possion function. There seems to 
be agreement also on the linear start of -P(cr) at small a, with a slope steeper that the WD 
function for the metallic regime. However the large s tail behavior is more controversial. 
Refs. [ |Tl|, |T6| claim that -P(o") has a Poissonian decay at large a 



P{a) ~ exp(-a(j), a » 1 , (6.5) 

with a ~ 1.9, wheras Refs. [ |T^, |T^] suggest a behavior in agreement with the plasma model 
result of Eq. (|6.4|) . We would like to emphasize that the numerical results published in Ref . [ 



T6| explicitly show good statistics for large values of a, and therefore we believe that they 
are reliable to extract the asymptotic behavior of -P(cr). 

The second important result provided by the numerical simulations is the existence of a 
linear term in the variance, as in Eq. ||6.3|| . These calculations do not exclude the presence 
of a power-law term of the kind ri{n)'^, which is however difficult to detect and quantify 
because of the presence of the dominant linear term.El. The coefficient of the linear term is 



shown to be ^ 0.27 in Ref. [ |T6l and 0.30 in Ref. [ ||]. 

The RME with weak confinement that we have considered in this paper is able to re- 
produce two of the main features seen in the numerical simulations of the critical statistics: 
the overall behavior of -P(cr)- with the linear rise at cr << 1 and the exponential decay at 
cr >> 1)- and the linear dependence of the number variance. In Fig. |1^ we plot the LSDF of 
the RME with logarithmic potential for A = 0A,f3 = 1, together with the critical LSDF of 



the Anderson model from Ref. [ 16]. The agreement between the two curves is spectacular in 
a very large energy range where -P(cr) varies by five orders of magnitude. Notice in the inset 
of the figure the behavior for large a of the tail of In [P(o")], which has apparently a linear 
slope. The parameter ~ 0.4 identifies, among all the possible members of the family of 
RME with logarithmic confinement, the ensemble which has the closest LSDF to the critical 
statistics. 

If we now compute the coefficient t] of the linear term in the variance for the RME with 
Ac ~ 0.4 [see Fig. figS], we obtain 0.32, which is consistent with the numerical results 
from the exact diagonalizations.lil'0 Thus the RME with logarithmic confinement is able 
to reproduce quantitatively the shape of the critical LSDF when A = A^ and at the same 
time provides an accurate estimate of the the leading order term of the number variance. In 
some sense, the "residual Possonian" properties of the critical statistics are well reproduced 
by this generalized RMT. The RMT does not provide the asymptotic power-law behavior 
( |6.2a| ) of the two-level correlation function found analytically, which, on the other hand, is 
also difficult to extract by direct numerical diagonalizations with reliable accuracy.^ Clearly 
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more work is necessary to determine if the good agreement shown here between the RME 
with weak confinement and the critical statistics is more than a furtuitous coincidence. It 
is however interesting and important that some of the properties of the correlations at the 
mobility edge can be correctly reproduced by such a simple ensemble. 



VII. CONCLUSIONS 

In this paper we have studied in detail the properties of families of random matrix 
ensembles (RME) that are invariant under similarity transformations but are characterized 
by a generalized level confinement. We have shown that the level statistics are affected 
by the confining potential when this is very soft. In particular, for a squared logarithmic 
potential, the statistical bulk properties are nonuniversal and deviate significantly from the 
Wigner-Dyson statistics of the Gaussian Ensembles, exhibiting a cross-over toward a more 
Poissonian behavior when an internal parameter is decreased. The U{N) invariant RME 
with logarithmic confining potential belongs, together the RME with a symmetry breaking 
term [ see Eq. ( |1.1| )], to a new universality class, distinct from Wigner-Dyson universality of 
classical RMT. 

We have shown that the nonuniversal behavior of the two-level correlation function for 
these RME can still be obtained within Dyson's mean field theory, generalized to the case 
weak confinement. 

We have performed Monte Carlo simulations to calculate several important statistical 
properties of the generalized RME that probe both short and long range correlations. 

The statistical properties of the RME with logarithmic confinement have strong similari- 
ties with the universal energy level statistics of disordered conductors at the metal-insulator 
transition. In particular, the probability distribution of the level spacings for a three di- 
mensional Anderson model at the critical point can be very well fitted, throughout a wide 
energy range, by the corresponding RME function for one particular choice of the internal 
parameter. Then for the same value of the parameter, this RME predicts a linear behavior 
for the level number variance, with a coefficient of proportionality close to the value obtained 
from numerical diagonalizations of the Anderson model at the critical point. 
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FIGURES 



FIG. 1. Density of states for the logarithmic potential for A = 0.5. The Monte Carlo (MC) 
results for (3 = 1,2 are plotted in small region around the origin e = with the MFT density, which 
corresponds to /3 = oo and diverges at e = 0. For /3 = 1, 2 and 4(not shown), the density is finite 
at e = 0. All the curve rapidly collapse on top of each other away from the origin. 



FIG. 2. The MC two-level correlation function R2{s' = 0,s) vs. s for the power-law potential 



with a = 0.2, [5 = 2. The solid line is the result of the CUE, R2{f') 
omitted the 6 function at the origin]. 
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FIG. 3. MC result for the two-level correlation function for the logarithmic confinement, 
Eq. ( |2.12| ), showing the existence of a "ghost" hole at s' = — s. The simulations are performed for 
[3 = 2 and A = 0.5 with = 101 particles. The reference particle is mobile around s ~ 24.4. The 
solid line in the inset corresponds to Eq. ([4.29| ). 



FIG. 4. MC results for the normal part of R2(n)is,s') [3 = 2, plotted together with the cor- 
responding exact expressions, Eq. (4.28), and with the GUE curve as a comparison. Notice the 
good agreement at small s' between numerical and exact results for A = 0.5 which becomes worse 
for A < 0.2. The fluctuations at s' > 1, more visible for A = 0.2, are due to finite-size fluctua- 
tions of the exact density, which the unfolding procedure via the MFT density cannot cure. These 
fluctuations are much smaller for [3=1. 



FIG. 5. MC result for the LSDF in the middle of the spectrum of the power-law potential with 
a = 0.2 and [3 = 2. For a ^ P{(j) vanishes like The exact GUE distribution is also 

plotted. 



FIG. 6. MC results for the bulk LSDF of the logarithmic potential with [3 = 1 and different 
values of A showing a cross-over betweem and the GOE and the Poisson distributions, also plotted. 
For A < 0.2 the LSDF stops approaching the Poisson function and tends to a single 5 function 
peak at o" ~ 1. Shown in the inset is the large a behavior of P{(t). Notice the logarithmic scale for 
the y-axis. The two dashed straight lines are fitting functions of the form exp{—a{A, [3)a) 



FIG. 7. The MC level number variance vs (n) for the power-law potential, for [3 = 1 and several 
values of a < 1. The dashed curve is the GOE result given in Eq. ( CT) . The energy windows 
containing {N) particles are centered at the origin e = 0. The nonuniversal "super- Wigner" 
behavior of the correlations at e = is responsible of the a-dependence of the constant term in the 
variance, which is otherwise equal to the that one of GOE. 
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FIG. 8. The relative variance var((?7-) / (n) vs {n) for the logarithmic potential for different values 
A and (3. The energy window are centered in the bulk of the spectrum and the variance is linear 
with (n). The constant straight lines are the slopes of the linear term of the exact solution of Ref. [ 
29| for (3 = 2 and corresponding A. 



FIG. 9. The MC variance vs. (n) for the logarithmic potential for the three symmetries and 
A = 0.5. The energy windows are now centered at the origin. Since the sum rule Eq. ( [4.32| ) is 
satisfied, the linear term in the variance (see Fig. |^ vanishes and the system becomes even more 
rigid than the GE. 



FIG. 10. LSDF of the RME with logarithmic confinement for (3 = 1, A = 0.4, plotted together 
with the energy-LSDF of the three dimensional Anderson model at the metal-insulator transition, 
taken from Ref. [ The inset shows the large a behavior in logarithmic scale for the y-axis. A 
straight line is the best fit for both curves. 
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